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ABSTRACT 

The melting of sea ice is correlated to increases in sea surface 
temperature and associated climatic changes. Therefore, it is 
important to investigate how rapidly sea ice floes melt. For 
this purpose, a new TempoSeg method for multitemporal seg- 
mentation of multiyear ice floes is proposed. The microwave 
radiometer is used to track the position of an ice floe. Then, 
a time series of MODIS images are created with the ice floe 
in the image center. A TempoSeg method is performed to seg- 
ment these images into two regions: Floe and Background. 
First, morphological feature extraction is applied. Then, the 
central image pixel is marked as Floe, and shape-constrained 
best merge region growing is performed. The resulting two- 
region map is post-filtered by applying morphological oper- 
ators. We have successfully tested our method on a set of 
MODIS images and estimated the area of a sea ice floe as a 
function of time. 

Index Terms — Classification, segmentation, region 
merging, multiyear sea ice floes, shape analysis. 

1. INTRODUCTION 

The sea ice cover has a significant effect on the Earth’s cli- 
mate system. It limits the exchange of heat, moisture, and 
momentum across the ocean- sea ice- snow-atmosphere inter- 
face. Therefore, sea ice is both an active participant of the 
Earth’s climate and a sensitive climate indicator [1]. It thus is 
very important to monitor sea ice evolution and develop meth- 
ods for automated analysis of satellite measurements. In this 
paper, our objective is to determine how rapidly a multiyear 
ice floe can melt. 

Most of the existing techniques for automatic ice floe 
tracking from Synthetic Aperture Radar or optical remote 
sensing imagery employ adaptive, or dynamic thresholding 
for distinguishing ice from water [2, 3]. In [4], Yu applied 
iterative region growing to the segmentation of sea ice im- 
agery from the RADARSAT-2 data. A watershed transform 
was performed, followed by labeling of the resulting homo- 
geneous regions and region merging to obtain an accurate 
segmentation map. 


In our study, two sensors onboard the NASA Aqua satel- 
lite have been used: Advanced Microwave Scanning Ra- 
diometer - Earth Observing System (AMSR-E) and Moderate- 
Resolution Imaging Spectroradiometer (MODIS). AMSR-E 
data was used to manually track the position of the ice floe 
of interest. Based on this information, a time series of the 
reprojected MODIS images was created with the ice floe in 
the central part of every image. These measurements made 
it possible to estimate both an area and a perimeter of a mul- 
tiyear Arctic sea ice floe as a function of time using a new 
multitemporal segmentation ( TempoSeg ) method. The 
proposed TempoSeg method is based on shape-constrained 
best merge region growing and proceeds as follows: 

First, a feature extraction by applying morphological op- 
erators is applied. Every pixel is considered as a separate 
region. The central image region is assigned to the class 
Floe , while all other regions are assigned to the class Back- 
ground. Second, hierarchical step-wise optimization segmen- 
tation [5] is performed, by iteratively merging two adjacent 
regions with the smallest Dissimilarity Criterion (DC). When 
any Background region is merged with the Floe region, the 
resulting region is assigned to the class Floe. The algorithm 
is converged for a given image when the image is segmented 
into two regions: Floe and Background. Finally, a morpho- 
logical post-filtering of the resulting map is applied. 

The important contribution of this work consists in includ- 
ing multitemporal, geometrical and shape analysis for com- 
puting a DC between two regions. Moreover, this method is 
applied on several images processed sequentially. 

The outline of the paper is as follows. In the next sec- 
tion, the data set used for experiments is described. Section 3 
presents the proposed TempoSeg segmentation method. Ex- 
perimental results are discussed in Section 4. Finally, conclu- 
sions are drawn in Section 5. 

2. DATASET 

In this study, we analyzed two months of AMSR-E and 
MODIS measurements. These sensors are onboard the Aqua 
satellite, which has a polar orbit enabling several measure- 
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Fig. 1. (a) Reprojected MODIS 1(227-1435) image, (b) Tem- 
poSeg segmentation map of 1(227-1435) (floe area = 24104 
pixels), (c) Reprojected MODIS 1(286-2045) image, (d) Seg- 
mentation map of 1(286-2045) (floe area = 19315 pixels). 

ments per day over the Earth’s polar regions. Our target 
of interest is a multiyear sea ice floe, which left the Arctic 
main sea ice pack in December 2007, continuously moved 
according to the Arctic ocean currents and melted in June 
2009. 

To track the position of the ice floe of interest during its 
lifetime, the AMSR-E data at 89 GHz mapped to a polar stere- 
ographic grid at 6.25 km spatial resolution were used. These 
data have been extracted from the “AMSR-E/Aqua Daily L3 
6.25 km 89 GHz Brightness Temperature Polar Grids” prod- 
uct distributed by the National Snow and Ice Data Center. 
In accordance with these measurements, a time series of the 
MODIS images was built with the ice floe of interest. 

In order to get the most accurate segmentation results, it is 
the best to use the highest MODIS spatial resolution band at 
250 m. Among the two high-resolution bands of the MODIS 
data, we used band 1 (0.620 - 0.670 pm) because it provided 
the best contrast between the sea ice and the ocean water. 
Level IB calibrated and geolocated MODIS swath data re- 
projected onto a polar stereographic grid was the main data 
source. The geolocation uncertainty is ^50 m [6]. 

All swath images did not provide opportunities to quantify 
the ice floe areas either due to extensive cloud cover, limited 
solar illumination, or weakness of contrast between the mul- 
tiyear ice floe and the surrounding young ice. We selected for 
our analysis M = 40 images with spatial dimensions of 500 
x 500 pixels, acquired on 18 different days during the two- 
months period from mid- August to mid-October 2008 (from 
the 227th day to the 286th day of 2008). Fig. l(a, c) and 2(a) 
show three images from the considered set. We denote each 
image by I m (D-T) = { i m (x,y ) G Z + ,x = [l..H],y = 
[1..W]}, m = [1..M], where D and T mean the day and the 


Fig. 2. Image 1(249-1715 ): (a) Reprojected MODIS im- 
age. (b) Segmentation map when no shape constraints and 
no post-filtering were applied (floe area = 23822 pixels), (c) 
TempoSeg segmentation map before post-filtering (floe area 
= 22903 pixels), (d) TempoSeg final segmentation map (floe 
area = 22853 pixels). 

time when it was recorded, H and W are height and width of 
each image, respectively (in our case, H = W = 500). For 
instance, the image I\ (227-1435) was acquired on the day 
227 at 14:35. 

3. PROPOSED METHOD 

The proposed TempoSeg segmentation method is based on 
shape-constrained best merge region growing procedure. It 
consists of the following steps: 

Input: At the input, M 1-band images I m (D — T) = 
{i m (x,y) e Z+,x = [1 -H],y = [l..W]},m = [1..M] in- 
eluding the Floe are given. The goal is to obtain for each 
image I rn a segmentation map S rn = { s m (x , y) G [0, 1], x = 
[ 1 . . H] , y m [ 1 . . W] } , where each pixel has label 1 if it belongs 
to the Floe , and 0 otherwise. The following parameters must 
be initialized: maximum number of pixels representing the 
floe Nmaxp, maximum height Hmaxp , maximum width 
Wmaxp of the flow, respectively, and its rectangularity at 
the previous time moment Rect° F . 

We define a height H eight (R) and a width Width(R) of 
the region R as the height and the width of the Minimum Area 
Rectangle (MAR) of an arbitrary orientation including R. A 
region rectangularity Rect(R) is computed as a ratio of the 
region area (i.e. number of pixels) and the area of the MAR 
enclosing R (i.e. Height(R) • Width(R)). 

Feature extraction: In order to enhance the difference 
between the intensity values of the Floe and the rest of the 
image, each image I m is filtered to remove small objects, as 
described in the following: 
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1. Morphological h- minima transform is applied to 7 m , 
which suppresses all minima in the image whose depth is less 
than a scalar h [7], resulting in an image 7^. 

2. 77- maxima transform is applied to 7^, which sup- 
presses all maxima whose height is less than h, resulting in 
an image 7^ [7]. 

3. Geodesic opening by reconstruction = 7 (7^) of 
the current image is obtained by reconstructing the image 
erosion e(Im) under 7^, using a disk with radius r as the 
structuring element (SE) [7]. This operation removes bright 
objects with a thinner support than the given SE, while pre- 
serving shapes of other objects. 

4. Geodesic closing by reconstruction 7^ = * s 

performed in the similar way, with a disk of radius r as a 
SE [7]. 

The obtained image 7^ is combined with the original 
image 7 m , resulting in a 2-band image I m = {i m (x,y) G 
Z 2 ,x=[l..H],y=[l..W}}. 

Hierarchical optimization: The segmentation by hierar- 
chical shape-constrained optimization of each image I m is se- 
quentially performed as follows: 

1 . m = 1 . 

2. Consider the image I m . Initialize the algorithm by 
assigning a new region label for each pixel. The one-pixel 
region situated in the center of the image [with coordinates 
(H/2+1, W/2+1)] gets a Floe class label L(R F ) = 1. All 
other regions get a class label L(Rj) = 0 , j = [1..W x 

3. Calculate the DC between all pairs of spatially adjacent 
regions { Rj,R k }, considering an eight-connectivity neigh- 
borhood and using the following algorithm (see Fig. 3): 


• Compute the square root of band sum mean squared 
error (BSMSE 1 / 2 ) criterion, which is based on mini- 
mizing the increase of mean squared error between the 
region mean image and the original image data. The 
BSMSE 1 / 2 between regions Rj and Rk with region 
mean vectors u j = ( Uj\^Uj^) T and = (u k i, u k 2 ) T 
and region sizes (number of pixels) rij and n k , respec- 
tively, is defined as 


DC(Rj,Rk) = BSMSE*{R j ,R k ) = 


TljTl k 

(rij + n k ) 


^(u jb u kk ) 


b= 1 


1/2 


( 1 ) 


• If either L(Rj) = 1 or L(R k ) = 1, i.e., one of the 
regions represents a floe, proceed as follows: 


a) If ( rij + n k ) > Nmaxp , DC(Rj , R k ) = 00 
(the upper maximum value of float). 

b) Otherwise, if (rij + n k ) > K • Nmaxp , 

0 < K < 1, perform shape analysis. The parame- 
ter K controls the scale at which a shape analysis 



Fig. 3. Flowchart of computing a DC between two regions 

Rj and R k . Rj k = Rj U R k , and / G {j, fc}, L(Rf ) = 1. 


is started. 

Let R jk = Rj U R k , and / G {j, fc}, L(i?/) = 1. 

• If (Height(Rj k ) > Hmaxp ) or 
(Width(Rj k ) > Wmaxp), then 
DC(Rj,R k ) = 00 . 

• Otherwise, if 

|i?ect(iijfc)— Rect™~ 1 \ < \Rect(Rf) — Rect™ -1 1, 
then DC(Rj,R k ) = 0.51X7(7^, i? fe ), 
else DC(Rj,R k ) = 2 DC(Rj,R k ). 

The idea behind this rectangularity condition is based 
on an assumption that the time interval between ac- 
quisitions of two consecutive images is typically short 
enough, so that the rectangularity values of the Floe re- 
gion within these two images are similar. 

4. Find the smallest DC value DC min. 

5. Merge a pair of spatially adjacent regions with DC = 
DC min- If one of the merged regions has a class label 1 , a 
new region gets a class label 1 , and its rectangularity value 
Rect(Rf) is recomputed. 

6 . If the image is segmented into two regions Floe and 
Background , go to step 7. Otherwise, update the DC values 
between the new region and all regions spatially adjacent to 
it, and return to step 4. 

7. Update parameters: 

• Red -p = 0.5 (Red™- 1 + Red™). 

• If m > 1: 
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Nmaxp = min [Nmaxp, 1-05 max (n™ 1 ,n 1 F )\. 

Hmaxp = min [Hmaxp, 

1.2 max(i7eig/i£(ii^ _1 ), i7eig/if(it”p ))]. 

W max p = min[VFmaxi?, 

8. m = m + 1. Go to step 2. 

Post-filtering: Each segmentation map S rn is post-filtered 
by applying the following morphological operators [7]: 

• Opening with a disk of radius r Q as a SE. 

• Area opening to remove all small Floe regions that may 
appear after morphological opening. 

• Closing using a disk with radius (r Q + 2) as a SE. 

4. EXPERIMENTAL RESULTS AND DISCUSSION 

We applied the proposed TempoSeg approach to the set of 40 
MODIS images described in Section 2. Based on the visual 
analysis of the first image in a time series I\ (227 — 1435) (see 
Fig. 1(a)), we initialized the parameters Nmaxp = 27000, 
Hmaxp = 250, W max p = 170 and Rect° F = 0.71 (com- 
pare to the floe characteristics we estimated by applying Tem- 
poSeg: n\ = 24104, Height (R],) = 232, Width^R)?) = 
145, Rect\ = 0.72). Other parameters were empirically set 
as h = 30, r = 15, K = 0.8 and r Q = 4. 

Fig. l(b, d) and 2(b-d) show some of the obtained segmen- 
tation images. From these images, it is possible to see that the 
proposed approach succeeded in accurately segmenting the 
ice floe of interest from the rest of the image. Fig. 2 demon- 
strates the improvement in segmentation results when includ- 
ing shape analysis into best merge region growing (compare 
Fig. 2(c) versus Fig. 2(b)) and when applying morphological 
post-filtering (compare Fig. 2(d) versus Fig. 2(c)). 

Fig. 4 depicts graphs of the Floe area and perimeter as a 
function of time. As expected, during the two-month period 
at the end of the summer season, both area and perimeter are 
decreasing over time. Some variability such as the local min- 
ima and maxima of the curves are a consequence of either a 
partial cloud cover over the Floe (minima) in some images, or 
weakness of contrast between the multiyear ice and the neigh- 
boring young ice floes (maxima). 

5. CONCLUSIONS 

We have designed a new TempoSeg method for multitempo- 
ral segmentation of multiyear sea ice floes from the MODIS 
data. The proposed technique is based on shape-constrained 
best merge region growing, and it segments each image from 
a time series into Floe and Background regions. We have 
applied this method to a set of MODIS images acquired in 
August-October 2008 and successfully estimated both an area 
and a perimeter of the floe of interest over the given time pe- 
riod. 




Fig. 4. (a) Ice floe area (number of pixels) as a function of 
time: Before post-filtering averaged for every day (blue); Af- 
ter post-filtering (green); After post-filtering averaged for ev- 
ery day (red), (b) Ice floe perimeter as a function of time, 
averaged for every day. 
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